Quantum ripples over a semi-classical shock 
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The evolution of an initially smooth spatial inhomogeneity in the density of a one-dimensional 
Fermi gas is well described by classical mechanics. The classical evolution leads to the formation 
of a shock wave: the density develops kinks in its coordinate dependence. We show that quantum 
corrections to the shock wave produce density ripples which run off the kinks. Despite their quantum 
origin, the amplitude and period of the ripples are expressed only in terms of classical objects derived 
from a smooth density profile. 
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In one-dimensional many-body physics, numerous 
long-wavelength properties of a quantum system are 
faithfully represented by the dynamical properties of a 
continuous liquid. Such representation is at the heart of 
a powerful bosonization method [IH3] . Small-amplitude 
perturbations of the liquid's density may be considered 
in a harmonic approximation. In that approximation, 
waves of density are linear, and there is no difference 
between the quantum and the classical dynamics of the 
density perturbations. At larger amplitudes, waves be- 
come nonlinear; it is not clear then, if the classical and 
quantum dynamics remain indistinguishable even in the 
long- wavelength limit. 

In a liquid made of fermions, higher density is associ- 
ated with higher Fermi momentum and higher particles 
velocities. Described classically, a higher-density soli- 
tary segment moves faster than its lower-density periph- 
ery. That leads to a wave overturn and the formation of 
shocks in the density p(x,t), i.e., points x 7 (t) where the 
x-dependence of p looses its analyticity. Formation of 
shocks can be seen within a classical continuous medium 
description in terms of the Riemann-Hopf equation • 
It is insensitive to the interaction strength, which can be 
set to zero (the free fermions case). The classical descrip- 
tion of propagation of shock waves in a gas of quantum 
particles is widely accepted across the fields of quantum 
physics from string theory [7j to cold atomic gases [8] . In 
case of free fermions, however, one may also investigate 
the fully-quantum evolution of the many-body wave func- 
tion, which is a Slater determinant of free-propagating 
single-particle states. 

We consider the shock formation in the semiclassical 
limit of the quantum evolution of free fermions. Our 
main finding is the appearance of oscillatory structure in 
dp/dx accompanying each of the shock points x^(t). This 
structure, quantum in origin, has a characteristic scale of 
variation fully determined by the solution of the classical 
problem. (In that respect, the free-particle many-body 
wave function associated with the shock formation bears 
some resemblance to the WKB wave function of a single 
particle in an external potential.) In classical physics, the 
formation of shocks [6] [9] is associated with the presence 
of velocity gradients. The inverse time to develop a shock 
can be estimated as ts ~ \dv/dx\, where v(x) is the initial 



velocity distribution. If the velocity gradient comes from 
a "bump" Ap in the density of fermions, then the cor- 
responding density gradient, dp/dx = (m/h)dv/dx be- 
comes large in the semiclassical limit, h — > at fixed 
dv/dx (hereinafter m is the mass of a fermion). Together 
with | dp/dx |, the total number of particles AN in the 
density bump scales as 1/h. We show that 1 /AN can be 
consistently used as a small parameter in the semiclassi- 
cal expansion. The characteristic length of the oscillatory 
structure around a shock point scales as (AiV)~ 2 / 3 . 

We are considering the time evolution of a solitary den- 
sity maximum p(x,t), which initially (at t = 0) has spa- 
tial extent Ax and amplitude Ap, 

p(x,0)=n + Apf(^-). (1) 

Here the dimensionless positive function /(s) reaches a 
maximum, /(0) ~ 1, at s = 0, has spatial extent |s| ~ 1, 
and decays outside that region, f(s — > ±oo) = 0, reveal- 
ing the background equilibrium density hq — Kp/n of a 
Fermi sea (\k\ < Kp) of spinless particles. We aim to 
describe the motion of semiclassical perturbations. That 
leads us to assume that many particles go into the cre- 
ation of the perturbation, AiV = ApAx 1. Our fur- 
ther assumption is that the density perturbation is small 
compared to the equilibrium density, Ap <C tiq . The two 
latter assumptions mean that KpAx AN ^> 1; there- 
fore, the density perturbation consists of particles and 
holes residing in the narrow vicinities Sk ~ KpAp/riQ of 
the Fermi points ±Kp. 

There are several characteristic time scales for the evo- 
lution of the density perturbation. First, it will de- 
compose into right- and left-moving modes, associated 
with particlc-hole excitations around each one of the 
Fermi points. The characteristic time associated with 
this process is t LR = (m/h)(Ax/n ). Next, there is 
a characteristic time over which the wave packet rep- 
resenting a single fermion spreads over a spatial scale 
Ax of the initial many-body state. This time is given by 
tq = (m/h)(Ax) 2 . On this time scale the evolution ceases 
to be semiclassical. Finally, in the course of the semi- 
classical evolution of the density perturbation, it may 
develop a shock, in analogy to the one-dimensional clas- 
sical hydrodynamics. This effect has been discussed, e.g., 
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in [TU]. The corresponding characteristic time scale is 
t s = (m/h)(Ax/Ap). 

Our assumptions AN 3> 1 and Ap <C n lead to the 
time scales hierarchy, t^ji -C ts <C Iq . Namely, first 
the perturbation splits into left and right moving modes, 
then classical shocks have a chance to appear in each 
of the modes, and at much later times the semiclassical 
nature of the evolution breaks down. We focus on the 
evolution occurring on time scale ~ ts- 

Since we are interested in the semiclassical regime, it 
is convenient to use the Wigner function, 



W{x,k,t) 



I* 



lky dy, 
(2) 



where |\&) is the initial (t = 0) state of the system. 

In equilibrium, |$) is a state where all single parti- 
cle levels with momenta between —KKp and hKp are 
occupied, while all other states are empty. Here Kp = 
irn . This state leads to a Wigner function given by 
W(x, k,t) = 6 {K F -k)d(k + K F ) . The Wigner func- 
tion is interpreted as a distribution in phase space of 
the Fermions, which, in equilibrium, occupy the region 
—Kp < k < Kp for all x, namely a band in phase space. 

We can think about a semiclassical state in which a 
single mode (either right- or left-moving) of density per- 
turbation is excited, as a state in which one of the Fermi 
wave numbers {—Kp or Kp) is allowed to vary with the 
space-time coordinates. From now on, and without loss 
of generality, we shall assume only a right moving mode, 
in which case the Fermi number Kp becomes time and 
space dependent, to be described by a function kp(x,t), 
while —Kp stays fixed. This heuristic picture translates 
into the following ansatz for the Wigner function: 



W(x, k, t) = 6(k F (x, t) - k)9(k + K F ). 



(3) 



We will see that Eq. (|3| is a good approximation to the 
Wigner function everywhere except the vicinities of the 
shock points. Knowledge of the Wigner function allows 
one to compute the density, as these two are related as 
follows: 



p(x,t) = / W(x,k,t) 



dk 
2tt' 



(4) 



Given initial conditions for the Wigner function, one may 
compute it at later times making use of the evolution 
equation, 

hk. 



d t W(x,k,t) = —d x W(x,k,t). 



(5) 



The proof of Eq. ^ consists of the application of the 
Schrodinger equation to the fermions in Eq. ^ and in- 
tegration by parts. The initial value problem for the 
evolution equation may be solved by the method of char- 
acteristics, which yields 



hk 

W(x,k,t) = W I x t,k,0 

m 



(6) 
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FIG. 1. Function kir(x,t) = mvF{x,t)/h at initial time 
(left), evolved according to Eq. |8} at an intermediate time 
(middle), and at a late time (right) where it becomes multi- 
valued. 



Applying Eq. (|6| to ansatz ^ one sees that kp(x, t) must 
satisfy: 



kp(x, t) = kp x kp(x, t)t, 

1 m 



(7) 



This equation determines kp{x,t) implicitly, given the 
initial conditions kp(x,0). 

The relation to the classical shock wave physics in one 
dimensional hydrodynamics may be realized by exclud- 
ing h from Eq. Q with the transformation vp{x,t) — 
(H/m)kp(x,t) and noticing that vp{x,t) is in fact a so- 
lution to 



dtVp(x, t) + vp(x, t)d x vp(x, t) = 0. 



(8) 



This Riemann (or Riemann-Hopf) equation [5] is well 
known to give rise to shock waves. 

In the following we make use of the inverse function, 
xp(k, t), defined as the solution to the following equation: 



kp(xF 1 1) 



(9) 



The evolution of kp(x,t) determines the evolution of 
Xp{x,t) as follows: 



(10) 



xp(k,t) — xp(k,0) H kt. 

m 

Therefore, the derivative dxp/dk — x' F (k,0) + ^t. If 
x^(£:,0) is negative, there is always a positive time at 
which dxp(k,t)/dk = 0. This equation, together with 
Eq. ( 10 ) determines the time evolution of points where 
\dkp(x,t) /dx\ = oo. This divergence is a manifesta- 
tion of the shock phenomenon; it is also termed as the 
'gradient catastrophe'. After the first time the function 
kp(x,t) displays an infinite derivative, the solution to 
Eq. ([7]) becomes multi-valued as shown in Fig. I] With- 
out loss of generality, we assume that the multi-valued 
region is bounded by two points, x~(t) and x+(t), be- 
tween which the function kp(x,t) has three branches, 
enumerated by kp\x,t) < kp\x,t) < kp\x,t). 

Regardless of the fact that kp(x,t) is multi- valued, 
initial conditions such as ([3| lead to a Wigner function 
which is 1 in a region of phase space between the curve 
kp(x,t) and the horizontal line k = —Kp, while the 
Wigner function is outside this region. Here the curve 
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FIG. 2. A typical density profile. The classical trailing and 
lead ing shock points are X—(t) and x+(t), respectively; cf. 
Eq. ( 16 1. The quantum ripple effect is amplified in the deriva- 
tive d x p(x,t). Inset: d x p/k 2 as a function of Sx — x — x+(t) 
in the region next to x + (t), see Eq. p4|. 



of kp(x, t) is understood to contain all the branches of the 
function, just as the curve described in the right panel of 

Fig. 00 

The integral in Eq. Q may easily be performed to 
obtain the density. The result is: 

2irp d (x,t) = Kp + kp\x,t) — kp\x,t)-\-kp\x,t), x € I 



2irp cl (x,t) = K F + k F (x,t), x£l, 



(11) 



where / = [x-(t),x+(t)] denotes the spatial segment of 
multi- valued behavior of kp . The superscript cl in p cl 
denotes that ( 11 ) gives only the result of classical theory, 



on which we want to improve in the vicinities of points 
x = x±. 

Next to X-(t) and x+(t) two branches of the function 

kp(x,t) meet. Explicitly, kp and kp 1 become equal to 
each other at x_(t): 



, F '(x-(t),t) = k { p\x^(t),t) = k F (t), (12) 
and similarly at x+{t): 

4 2) (x+(t),<) - 4 3) (^+ (*),*) = 4w ■ (is) 

The functions kp (t) are solutions of the equation 

f 1 2") 

dxp{k,t)/dk = 0. The branches k F ' display a square- 
root behavior next to X-{t): 



k { p 1] (x,t) - kp{t) = ±a~{t)y/x-x-(t) + 



(14) 



(3) 

while k F is regular next to this point. At the same time 
next to x + (t) we have: 

kf 2) {x,t)~ k F (i) = ±a+ (t) y/x + (t)-x + ..., (15) 

while now kp^ (x, t) is regular. The ellipsis denote regular 
terms. 



Near the branch meeting points, p{x, t) takes the ap- 
proximate form: 



p Cl (x,t) - p Cl ( X± (t),t) ~ ±^My/ T { X - X± ( t )) . (16) 

7T 

These square-root kinks (which may be observed in 
Fig. [5] next to X- (t) and x+ (t) , albeit rounded by quan- 
tum corrections) point to the non-analytical behavior 
of p(x,t) at x = x±(t). The corresponding coeffi- 
cients are determined by initial conditions, [a (t)] 2 = 
2{df,xp)~ 1 . According to Eq. (10), the derivative d%Xp = 
d 2 xp{k, 0)/dk 2 at fixed k is time- independent. To find 
a (t) one should evaluate the said derivative at k = 
k F (t). 

Our goal is to evaluate quantum corrections to the clas- 
sical density profile (16). To that end we start with in- 



troducing a family of many-body coherent states of free 
fermions |llj with smooth density profiles. We are inter- 
ested in times t 3> thRi so it is sufficient to include only, 
say, right-movers in the consideration. The correspond- 
ing fermion density operator is p R = (x)ip^ R \x) , 
and the field operators contain only Fourier harmonics 
with k > 0, 



i> {R \x) = / *P k e 

Jk>0 



ikx dk 
2tt 



The coherent states 

|vf/) = e l f^( x )p R ( x ) dx \Q^ 



(17) 



(18) 



are parametrized by the function <5>(x) with a transparent 
meaning, 2irp cl = Kp + 2ir$'(x). 

The standard bosonization methods \TT\ yield 

p i(*(xi)-*(a5 2 )) 

(*\^(xi)iP(x 2 )\V) = — -, (19) 

i2tt(xi — X2) 

where we customarily dispensed with the additive <&- 
independent term which formally vanishes in the large- 
Kp limit, Kp\x\ — x%\ — > 00. 

To proceed with computing the Wigner function 
W(x, k,0) we perform a gradient expansion in the ex- 
ponent of Eq. (19), to obtain: 



W t ( x+ | j0 )^(x-|,o)|*) 



I \^y 3 +k F (x)y\ 



i2-Ky 



(20) 



The Fourier transform with respect to y yields an ap- 
proximation to the Wigner function [12 : 

W(x, k, 0) fa Aix (2 2/3 k(x, 0)(fc - k F (x, 0))) , (21) 

where 



Ai 1 (x) = / Ai(s')^'= / e l (^ + ^ dq 



2iriq ' 



(22) 
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and k is the t = value of 

K(x,t) = \d 2 x k F {x,t)/2\- 1 ' 3 . 



(23) 



The function Aii(s) changes monotonically from 1 at s — ¥ 
— oo to at s — > oo on the scale \s\ ~ 1. Therefore at 



large n(x, 0) the Wigner function of Eq. (21 ) approaches 
its classical form for the right-movers, 



W(x,k,0) « 9{k F (x) - k). 



This confirms the ansatz pj), where the second factor 
9(k + k F ) is now missing, as accounting exclusively for 
the right-movers formally corresponds to sending the left 
Fermi point, —K F , to — oo. 

Estimating n(x, 0) with the help of Eq. ([lj as n ~ 
A J o _1//3 Acc 2 / 3 , we see that smearing of the Fermi step- 



function, as described by Eq. (21 1, is smaller by the pa- 
rameter (AiV)~ 2 / 3 than the shift ~ AN/ Ax of the Fermi 
wave vector k F (x) from its equilibrium value. The same 
semiclassical parameter is needed to justify the gradient 
expansion which we used to derive Eq. (21 1. To assess 



the accuracy of that approximation, we may consider the 



correction to Eq. (21 ) resulting from the sub- leading term 



in the gradient expansion. The correction has the form: 



8W 



K 5 k F (x) 
2 2 / 3 5! 



Ai 



i"" (2 2 / 3 n{x,0){k - k F (x,0))\ . (25) 



Now we again use Eq. (Tlj), to find k^ 



ApAx- 



addition to the est ima te for n(x, 0). Substitution of these 
estimates in Eq. ^ yields 5W/W ~ AN~ 2 / 3 at \k - 
k F (x)\<l/n(x,0). 



We wish now to propagate the initial conditions ( 21 ) in 
time. It will be more convenient, however, first to switch 
to a representation involving only the inverse function, 
x F (k,t), defined in ([9]), rather than k F (x,t). We first 
focus on the factor k— k F (x : 0) in |2l"] ), aiming to get rid 
of the function k F and replace it by x F . Making use of 
([9]) to replace k by k F (x F (k, 0), 0) and Taylor expanding 
the function k F around x F {x, 0) we find: 



k F (x,0) 



x — x F (k, 0) 
x' F (k,0) 



(26) 



where we have also used k' F (x F (k, 0), 0) = l/x' F (k, 0) 



We substitute the first term of expansion (26 1 into (21 1. 
This is justified at x — x F (k, 0) <^ Ax, and at the same 
ti me su fficient to allow the argument of the Aii function 
in (21 ) to vary in the range \s\ ~ (AN) 2 / 3 > 1. In that 



range, the Wigner function closely approaches its limits 
1 and at negative and positive values of the argument, 
respectively. 

We turn our focus now to k(x,0) appearing in (21 1. 
Since n, as defined in (23), contains d 2 k F (x, 0), we need 



to replace it by an expression involving only x F (k, 0). To 
this end we use the following classical differential identity: 



d 2 x k F (xM 



■1 , 



=„(k) 



-;(fc,o) 

x' F (k,0) 3 ' 



(27) 



Combining Eqs. (26) and (27) we may re-write Eq. (21) 
as 



W(x, k, 0) « Aii (2 2/3 R{k)(x - x F {k, 0)) 



where 



k{k) = h(k,t) = \dlx F (k 1 t)/2\ 



-1/3 



(28) 



(29) 



(24) is independent of time, according to Eq. ( 10 ) 



We are now ready to compute the Wigner function 
W(x, k, t) and quantum corrections to density at later 
times. Combining Eqs. (10), and (28) we easily find: 



W{x, k, t) w Aii (2 2/3 k(k)(x - x F {k, t)) 



(30) 



Given the validity of Eq. ( |30| for all times, it remains 
only to integrate it over k to obtain the density. Using 
the representation (22) of Aii, one arrives at: 



d x p{x,t) 



-q 3 +q{x-x F (k,t)) 



]dq 



dk 



2ir 2ir 



(31) 



Next we implement a gradient expansion of function 
x F (k,t) around a branch meeting point, x±(t). Using 
this expansion in Eq. (31) and changing variables there 
to K = k + I and Q = k — | , we derive 

d x Sp ± (x,t) 



(32) 



l»f (*J ffljl ,,, , ! 



) + (x-x ± (t))(K-Q)\ dQ dK 

2tt 2tT 



Note that here we are computing only the singular con- 
tribution dp (x,t) next to the branch meeting points, 



8p ± (x,t) = ±(p(x,t) - p cl (x ± (t),t)) 



(33) 



rather than full density p(x,t). Integration over Q and 
P in Eq. (32) factorize, leading to the result: 



d xP ± (x,t) » ^(t) 2 [Ai^^fai)] 1 



(34) 



Here the distances, 8x± — =p(x — x±(t)), are measured 
from the classical shock points, and inverse length scales 
k (t) are related to the parameters of the classical so- 
lution K±(t) = (a ± (t)) 2 / 3 , cf. Eq. ((l6j). The density 
gradient (34) exhibits pronounced oscillations, "quantum 



ripples" in the vicinity of the classical shock points, see 
the inset in Fig. [2j At times t >t$, the distance between 
the shock points is x+(t) — X-(t) ~ Ax. The characteris- 
tic length scale for the ripples 1/K ± (i) ~ Ax /(AN) 2 / 3 
is parametrically smaller. Making use of the identity 
Ai'(x) = xAi(x), one may integrate Eq. (34) to find 



6p ± (x,t) (35) 
~ (t) { [Ai'(K± {t)6x ± )] 2 - («± (t)Sx ± ) Ai 2 (k± (t)Sx± ) } . 



Equation (J35J) together with the derivative (34) is the 
main result of this Letter. A typical plot of p(x) is shown 
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in Fig. [5J with the inset showing the graph of d x 8p ± (x) 
around a branch meeting point. 

To conclude, we found quantum corrections to the clas- 
sical shocks forming in the course of evolution of a one- 
dimensional Fermi gas. The leading quantum effect is 
the formation of density ripples in the vicinity of the 
shock points. The quasiclassical parameter controlling 
the spatial structure and strength of the ripples scales as 
(AN)- 2 / 3 with the excess number of particles AN in the 
evolving nonlinear wave. 
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